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Abstract 



Statistical physics is employed to evaluate the performance of error- 
correcting codes in the case of finite message length for an ensemble of 
Gallager's error correcting codes. We follow Gallager's approach of upper- 
bounding the average decoding error rate, but invoke the replica method to 
reproduce the tightest general bound to date, and to improve on the most 
accurate zero-error noise level threshold reported in the literature. The rela- 
tion between the methods used and those presented in the information theory 
literature are explored. 
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Many of the problems addressed in the Information Theory (IT) literature show great 
similarity to those treated in statistical physics. One of the main areas where these links 
are particularly strong is that of digital communication and coding theory; these links have 
been recently examined in the area of Low Density Parity Check (LPDC) and turbo || 

error-correcting codes. It is only natural to expect that some relations between the ana- 
lytical methods used in the two disciplines will emerge, and that advances in one could be 
employed to improve results in the other. In this Letter we focus on such an example. We 
utilize the replica method of statistical physics to assess the performance of Gallager's error 
correcting code in the case of finite message length, generalizing an established method in 
the IT community. The analysis reproduces the tightest general bound to date, but more 
importantly, it provides exact results to specific code constructions. 

Error correcting codes play a vital role in facilitating reliable data transmission, ranging 
from cellular communication to data storage on magnetic media. In a general scenario, the 
N dimensional Boolean message £ e {0, 1}^ is encoded to the M(> N) dimensional Boolean 
vector Zo, and transmitted via a noisy channel, which is taken here to be a Binary Symmetric 
Channel (BSC) characterized by flip probability p per bit; other transmission channels may 
also be examined within a similar framework. At the other end of the channel, the corrupted 
codeword is decoded utilizing the structured codeword redundancy. 

The block error rate Pe, defined as the probability for a decoding error, serves as a 



performance measure for the success of the coding method. In his seminal work |T3| , Shannon 
showed that the error rate can vanish for code rates R below the channel capacity in the 
limit N, M — > oo; in the case of the BSC and unbiased messages R = N/M < 1 — H 2 (p), 
where H 2 {p) = —p\og 2 p— (1—p) log 2 (l —p)- The upper bound, for infinitely long messages, 
is often termed Shannon's limit to the error correcting ability. Evaluating Pg for practical 
codes of finite length became one of central topics in IT. 

For maximum likelihood (ML) decoding where the most probable message given the 
possibly corrupted codeword defines the message estimate, it is believed that Pe of the best 
code scales as exp[— ME(R)\. The non-negative exponent E(R) is termed reliability function 



(RF); it becomes positive below the channel capacity defining the sensitivity of the optimal 
error rate to the message length, complementing Shannon's result. 

Unfortunately, assessing the RF directly is generally difficult. Instead, Gallager's pow- 
erful method || bounds E(R) from the below utilizing the inequality 



which holds for any arbitrary ML estimation, inferring a binary vector x after observing a 
vector y, and a positive variable p>0. 

The average error rate Pe for a certain ensemble of codes is greater than the ensemble 
minimum. Therefore, averaging the RHS of Eq.(|]) over the ensemble, one obtains an upper- 
bound to the minimum error rate that scales exponentially with M for large but finite N 
and M, exp[—ME av (p,R)]; the exponent E av (p,R) serves as a lower-bound of E(R). One 
can tighten the lower bound by maximizing E av (p, R) with respect to p>0. 

Evaluating E av (p, R) is also difficult (except for pGlV). The strategy used by Gallager || 
is to further upper-bound the RHS of Eq.([TJ) utilizing Jensen's inequality (x p ) < (x) p , which 
holds for any < p < 1 with respect to the expectation over any arbitrary distribution of a 
positive variable. The added inequality presumably makes the bound looser. It is therefore 
surprising that maximizing the exponent with respect to p G [0, 1] in the ensemble of all 
random codes having the same rate R, which results in the random coding exponent E r (R), 
provides an exact evaluation of the RF for high R values. 

However, the bound by E r (R) becomes loose once the optimal value of p reaches the 
upper limit of the interval, i.e., p = 1 (corresponding to Bhattacharyya's bound). It is 
not clear whether Jensen's inequality or Gallager's inequality ([I]) is responsible for this 
breakdown. Moreover, it is unclear how to devise a similar method for deriving bounds for 
other (non-random) codes, a question of high practical significance. 

In this Letter we demonstrate how the methods of statistical physics may be employed 
to obtain tighter bounds for specific codes. This is carried out by a direct evaluation of 
E av (p,R) for the ensemble of Gallager error-correcting codes ||. This (linear) code was 
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rediscovered only recently [J/J, showing outstanding performance, competitive to other state- 
of-the-art techniques. It is characterized by a randomly generated (M — N) x M Boolean 
sparse parity check matrix H, composed of K and C (> 3) non-zero (unit) elements per 
row and column, respectively. Encoding the message vector is carried out using the 
M x N generating matrix G T , satisfying the condition HG T = 0, where z = G T £ (mod 2). 
The M bit codeword z$ is transmitted via a noisy channel, BSC in the current analysis; the 
corrupted vector z = z + £ (mod 2) is received at the other end, where £g {0, 1} M represents 
a noise vector with an independent probability p per bit of having a value 1. Decoding is 
carried out by multiplying z by the parity check matrix H, to obtain the syndrome vector 
J = Hz = H(G T £ + £) =HC, (mod 2), and to find the most probable solution to the parity 
check equation Hn = J (mod 2) , for estimating the true noise vector £. One retrieves the 
original message using the equation G T S = z — n (mod 2); S to estimate of the original 
message. 

To facilitate the analysis we map the Boolean (0, 1) variables onto the binary (±1) 
representation. The binary vectors n and J, represent the noise estimate and syndrome 
vectors respectively; the latter is generated by taking products of the relevant noise bits 
Jft = Cii M --Cix M ; where the indices i^, ..^k^ correspond to the nonzero elements in row // of 
the parity check matrix H. 

The similarity between error- correcting codes and physical systems was first pointed out 
by Sourlas mapping a simple Boolean code onto Ising spin models with multi-spin 
interactions. We recently extended his work to more practical parity check codes 0. We 
employ a similar formulation using the Hamiltonian 



M 



H(n;J)= 7 ^ J D g 5 \j g - - ftm - F^n, , (2) 
g \ ieg J i=i 

to evaluate the joint probability for J and n 

Here, Q = (z'i, .., ik) runs over all combinations of K indices out of M; Jg = Yiieg an d the 
sparse tensor Dg becomes non-zero (unit) only when all indices in Q correspond to non-zero 



(unit) elements in a certain row of the parity check matrix H . Taking 7 — > 00 enforces the 
parity check equation. The additive field F— (1/2) In [(l—p)/p\ corresponds to the true prior 
probability in the Bayesian framework, reflecting the flip rate p. The inverse temperature 
(3 is introduced to emphasize the link with the statistical mechanics formulation and is 
generally fixed to (3 = 1 unless specified otherwise. 

One can then use (f|) to evaluate Pe from ([]]) by calculating the bound without invoking 
Jensen's inequality. The first part of the Hamiltonian @ is invariant under gauge trans- 
formations of the form rij — > riiQ, and Jg — > Jg Tiieg d — 1> which decouple the correlation 
between the dynamical vector n and the true noise £. Rewriting the Hamiltonian one ob- 
tains a similar expression to Eq. (fj) apart from the last term on the right which become 



Quenched averages over the ensemble of codes is carried out with respect to the current 
random selection of the sparse tensor D and the noise vector, which eventually results in 
a similar procedure to the replica method in statistical mechanics. This gives rise to a set 
of order parameters q a ,fi,...^ — jf J2i=i %i n f n i--- n l > where a, (3 . . . represent replica indices, 
and the variable Z% comes from enforcing the restriction of C and L connections per index 
respectively as in |J . This interesting similarity between Gallager's method and the replica 
method has been pointed out by Iba in [|J. 

To proceed further one has to make an assumption about the order parameter symmetry. 
As a first approximation we assume replica symmetry (RS) in the following order parameters 
and the related conjugate variables 



where I is the number of replica indices, q and q are normalization variables (7r(x) and tt(x) 
are probability distributions). Unspecified integrals are over the range [—!,+!]. 



shown that in the limit of large M this becomes identical to the full summation in the 
non-ferromagnetic phase, where 7r(x) 7^ 5(x — 1) and tt(x) 7^ S(x — 1). Then, one obtains 




(4) 



Originally, the summation Tr| n _^(-) excludes the case of n 7^ however, it can be 
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the expression 



E av {p,R) = - — In 



Tr P—(JX)\ Tr P~p(J,n 

V {J,0 Vm^C} 

F \\ 
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■ In (2 cosh F) — In I 2 cosh 
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in z£ F c,£>; 



(5) 



S>li+p>- 

where Z NF (£, L>; ^-^j) denotes the partition function Tr n lim-^oo exp[— (37i] in the non- 
ferromagnetic phase for a system with an effective additive field F/(l + p). Averages 
(■}^I_f_ d are over the distribution P(C; j^) = expf^-^ X^i d}/ (2 cos h(i+^)) an d the 

uniform distribution of D. Extremizing (Z^ F (CiD; jf^)),- F with respect to the order 

^ > I i+p >d 

parameters q, q, tt(-) and 7f(-)> under the replica symmetry ansatz (^), one obtains for the 
final term in (|5f) 



: Ext* ( ^ / n f 1+u ^ Xt 
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(6) 
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where Ext* denotes extremization which excludes the ferro- magnetic solution and (-)^ f is 



i+p 



over P(C, t^)- 

Before proceeding any further, we would like to mention some general properties of 
E av (p,R). From Eqs. (||) and (|6]), it can be shown that \im p _> E av (p, R) = and 
d 2 E av (p, R)/dp 2 < 0. This implies that Max p>0 E av (p, R), becomes positive if and only 
if dE av (p, R)/dp\ p=0 > 0, for which lirriM Pe = holds. Therefore, the zero error thresh- 
old, defined as the critical flip rate below which the average error rate vanishes as M — * oo, 
is obtained by the condition dE av (p, R)/dp = 0. From (||), this becomes 



F tanh F - — (In Z NF (C,D; F)) C]FD = 0. 



(7) 



The second term is the averaged free energy for the Hamiltonian (0) with respect to the 
quenched randomness £ and D, in the non-ferromagnetic phase. Employing the ferromag- 



netic gauge |T(| one obtains the following expression for the ferromagnetic free energy (where 



Pe — O)' (1/M) (InZp (£, D; F))^ FD = FtanhF. Since the correct prior information about 
the flip rate p is used in the calculation, these two free energies are actually obtained in Nishi- 
mori's finite decoding temperature ((3=1) ||I2],II,|I0,|5|| for which the bit error probability is 
minimized. By satisfying (|7|), the zero error threshold for ML decoding, which corresponds 
to the zero temperature limit ((3 — >oo) [0,0, is determined by the phase boundary between 
the ferromagnetic and non-ferromagnetic phases at f3=l. 

Using the ferromagnetic gauge provides insight into the physical properties of the system. 
As the internal energy per bit in the non-ferromagnetic system is — FtanhF under Nishi- 
mori's condition, Eq. (0) implies that the entropy of the non-ferromagnetic phase vanishes 
at the phase boundary for (3 = 1, suggesting that this phase exhibits a replica symmetry 
breaking (RSB) at lower temperatures in general, and at (3— >oo in particular. In this sense, 
the zero-error threshold prediction obtained from Gallager's method and ML decoding, is 
surprising as it provides information about the ferro/non-ferro phase boundary at (3 — > oo 
which is not easily obtained via the methods of statistical physics due to RSB effects. This 
argument can be extended to the case of general (3 > 1, as will be presented elsewhere. 

An analytical expression to E av (p, R) can be obtained in the limit K,C — > oo, keeping the 
code rate R = 1—C / K finite; for the non- ferromagnetic solution one then obtains q = 2 p l K , q = 
2 p(i-i/a-) ) 7r ( a .) = < j( x ) anc i 7f(x) = (l/2)(l+tanhF)5(x-tanhF)+(l/2)(l-tanhF)(5(x+tanhF). 
Using Eqs. (H) and (||), one obtains the explicit expression E av (p,R) = In 2 cosh F — (1 + 
p) In (2 cosh +p(l — R) In 2. In addition, there exists another solution for p > 1, q = 
2 i/* q = 2 i-i/a 7r(^c) = (l/2)5(x - l)+(l/2)5(x+l) and 7? (x) = (l/2)<J(x-l)+(l/2)<J(x+l) 
providing E av (p, R) = In 2 coshF — In (2 coshF+2 cosh (j^F^ + (1 — i?)ln2. Employing a 
method similar to that in 0H, it can be shown that both RS solutions are locally stable 
against perturbations to the replica symmetric solution. 

The relation between E av (p,R) and the entropy of non-ferromagnetic solutions 5nf 



dE av (p } R) ( Z nf(C,At^) 5 nf(CAt^)>C 



1 i+p ' 



suggests another type of RSB, indicated by the negative entropy. This implies that the 



Max E av {p,R)-- 

p>0 



entropy of the non-ferromagnetic RS solutions vanishes at p = p*{R) which maximizes 
E av (p, R); and the tightest lower bound of E(R) is therefore obtained at the RSB transition, 
which can be calculated from the locally stable RS solutions. 
Solving the maximization problem one obtains 

In 2 cosh F- (1-/2) In 2 F>2F*(R) 
-ln(2coshF+2) , 

ln2coshF-(l-i?)ln2 2F*(R)>F>F*(R) (8) 
-Ft&ahF*(R) , 

, otherwise 

where F*(R) is the solution of the equation In 2 coshF* — F* tanhF* — (1 — R) In 2 = 0. 
The position of the maximum is given as p*{R) = 1 for F > 2F*(R), F/F*(R) — 1 for 
2F*(R) >F> F*(R) and 0, otherwise. Using the relation between F and p, this indicates 
that E(R) becomes positive if and only if R< l — H^ij))-, which corresponds to Shannon's 
limit. 

Equation @ is identical to the random coding exponent E r (R) obtained in the IT liter- 
ature [Q, although one should emphasize the main differences between the two approaches: 
a) Strating from Gallager's inequality ([]]) we directly average over the ensemble while the 
E r (R) result is obtained by invoking Jensen's inequality, b) Our result is obtained for an 
ensemble of a specific code. 

With some hindsight, this is not very surprising as Gallager codes become similar to 
random codes in the limit K, C — > oo |7j|| ; this also implies that using Jensen's inequality 
does not produce a looser bound as initially thought. 

To get a tighter bound for low R values we employ a refined inequality, upper-bounding 
the ensemble minimum of Pe by ^Tr^j ^ P^+p(J, £) (Tr^j n^} P~(J> n )) ^ ^ (p > 
0, m > 0), as in ([!]). A similar calculation along the lines described here (details will be 
shown elsewhere) provides the expurgated exponent bound 0] result for low R values (see 
Fig.l); this links our results to the best bounds reported in the IT litereture to date. 
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Without trivializing the results obtained in the case of K, C — > oo, the main achievement 
of our approach is the ability to investigate analytically the performance of Gallager (or 
similar) codes of finite K and C. To demonstrate the accuracy of the bounds obtained we 
examine the case of K = 6 and C—3. We numerically evaluated E av (p, R) (^) for p— 0.0915, 
a recent highly accurate estimate of the error threshold for this parameter ]]J, and for 
p = 0.0990, which is the threshold predicted by our analysis. The numerical results were 
obtained by approximating tt(-) and tt(-) using 10 6 dimensional vectors and iterating the 
saddle point equations until convergence. The results are shown in the inset; they indicate 
that Max p >o E av (p, R) ~ 1.0 x 10~ 4 > for p — 0.0915 while E av (p,R) is maximized (to zero) 
in the vicinity of p = for p = 0.0990, suggesting a tighter estimate for the error threshold 
than those reported in the IT literature. 

In summary, we have developed a method to tightly upper-bound the dependence of 
the decoding error rate on the message length for Gallager codes. In the limit of infinite 
connectivity our result collapses onto the best general random coding exponents reported in 
the IT literatures, the random coding exponent and the expurgated exponent for high and low 
R values respectively. The method provides one of the only tools available for examining 
codes of finite connectivity; and predicts the tightest estimate of the zero error noise level 
threshold to date for Gallager codes. It can be easily extended to investigate other linear 
codes of a similar type and is clearly of high practical significance. 

We demonstrated how the methods of statistical physics may complement and improve 
results obtained in the IT literature. These methods are applicable to a broad range of 
problems, especially within the sub-field of coding, and may be instrumental in improving 
existing results; some of these studies are already under way. 
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FIGURES 




FIG. 1. Lower-bounds on the reliability exponent E(R) obtained for p = 0.01 in the limit 
K, C — > oo. Our method produces the same result as the random coding exponent E r {R) (solid 
line) which provides an excellent bound for R> Rj,. For low R<R a values the bound becomes loose, 
and a better result (dashed line), identical to the expurgated exponent bound, is obtained (see text) 
by employing a refined inequality in (|l|). Inset - The exponent E av (p,R) obtained numerically for 
a choice of finite parameters K = 6 and C = 3 (R = l/2). Symbols and and standard deviations are 
computed using 50 numerical solutions. Curves are obtained via a quadratic fit. For p = 0.0915, 
p*(R) — 0.02, suggesting that this flip rate is still below the threshold. Left of the peak, the RS 
solution (thin broken curve) is unstable. For p = 0.0990, our predicted threshold, the maximum 
E av (p,R)~0 is obtained at p — 0, implying that this is the correct threshold. 
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